function [quality, r] = quality_measure(method, varargin)
% varargin = im1, im2, scale
% varargin = {im1, im2}, scale
% varargin = {im1}, {im2}, scale

nIn = length(varargin);
if length(varargin{end}) == 2
    nIn = nIn - 2;
    [scale, shift] = varargin{end-1 : end};
else
    nIn = nIn - 1;
    shift = [0 0];
    scale = varargin{end};
end

switch lower(method)
    case 'psnr'
        hQuality = @(x, y) calc_psnr(x, y, scale, shift);
    case 'ssim'
        hQuality = @(x, y) calc_ssim(x, y, scale, shift);
end

switch nIn
    case 1
        [quality, r] = calc_quality(varargin{1}(:, 1), varargin{1}(:, 2), hQuality);
    case 2
        [quality, r] = calc_quality(varargin{1:2}, hQuality);
end

%%
function [psnr, r] = calc_psnr(im1, im2, scale, shift)
im1 = srimg.preProcess(im1, scale, shift(1));
im2 = srimg.preProcess(im2, scale, shift(2));

r = double(im1) - double(im2);
rmse = sqrt(mean(r(:).^2));
psnr = 20 * log10(255 / rmse);

%%
function [ssim, r] = calc_ssim(im1, im2, scale, shift)
r = [];
im1 = srimg.preProcess(im1, scale, shift(1));
im2 = srimg.preProcess(im2, scale, shift(2));

% ssim_val = ssim(im1, im2);
ssim = srimg.ssim_index(im1, im2);

%%
function [quality, r] = calc_quality(x, y, hQuality)
% if iscell(x) && iscell(y)
%     [quality, r] = @(c1, c2) cell2mat(cellfun(hQuality, c1, c2, 'UniformOutput', false));
%     return
% end

dim = size(x, 3);
if dim == 3
    [quality, r] = hQuality(x, y);
    return
end

quality = zeros(dim, 1); r = cell(1, dim);
for i = 1 : dim
    [quality(i), r{i}] = hQuality(x(:, :, i), y(:, :, i));
end
